######## Attitudes to crime and punishment in England and Wales, 1964-2023
######## 002. Analysis_Demographics

#### Set up

library(tidyverse)
library(gt)

Data <- readRDS(file = "/Users/matteo/Desktop/Prison Polling/AttitudesToCrimeAndPunishment.RData")

#### Demographic correlations

analyse_data <- function(data, demog, var, smoothing = TRUE) {
  data %>% filter(Demographic == demog) -> data
  dyads <- DyadRatios::extract(data$Varname, data$Date, data$Index, 
                               smoothing = smoothing)
  return(tibble(Trend = dyads$latent1,
                Year = dyads$period,
                Topic = var,
                Demographic = demog))
}
corr.results <- matrix(nrow = 9, ncol = 1)

Data[[1]] %>%
  mutate(Date = as.Date(Date)) %>%
  filter(Date > as.Date("1980-01-01")) %>%
  select(Varname, Date, Index, Poll, Demographic, n) %>%
  mutate(Demographic = ifelse(Demographic %in% c("16-24","25-49"),"Under 50",
                              ifelse(Demographic %in% c("50-64","65+"),"Over 50",
                                     Demographic))) %>%
  group_by(Date, Demographic, Varname) %>%
  summarise(Index = weighted.mean(Index, n)) -> temp

results <- analyse_data(data = temp, demog = "Over 50", var = "FearCrime", smoothing = FALSE)
analyse_data(data = temp, demog = "Under 50", var = "FearCrime", smoothing = FALSE) %>%
  bind_rows(results, .) -> results
analyse_data(data = temp, demog = "Male", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp, demog = "Female", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp[temp$Date >= as.Date('1991-01-01'),], demog = "White", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp[temp$Date >= as.Date('1991-01-01'),], demog = "Non-white", var = "FearCrime") %>%
  bind_rows(results, .) -> results

corr.results[1,] <- cor.test(results$Trend[results$Demographic == "Over 50"], results$Trend[results$Demographic == "Under 50"])$estimate
corr.results[2,] <- cor.test(results$Trend[results$Demographic == "Male"], results$Trend[results$Demographic == "Female"])$estimate
corr.results[3,] <- cor.test(results$Trend[results$Demographic == "White"], results$Trend[results$Demographic == "Non-white"])$estimate


Data[[2]] %>%
  mutate(Date = as.Date(Date)) %>%
  filter(Date > as.Date("1980-01-01")) %>%
  select(Varname, Date, Index, Poll, Demographic, n) %>%
  mutate(Demographic = ifelse(Demographic %in% c("16-24","25-49"),"Under 50",
                              ifelse(Demographic %in% c("50-64","65+"),"Over 50",
                                     Demographic))) %>%
  group_by(Date, Demographic, Varname) %>%
  summarise(Index = weighted.mean(Index, n)) -> temp

results <- analyse_data(data = temp, demog = "Over 50", var = "FearCrime")
analyse_data(data = temp, demog = "Under 50", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp, demog = "Male", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp, demog = "Female", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp[temp$Date >= as.Date('1991-01-01'),], demog = "White", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp[temp$Date >= as.Date('1991-01-01'),], demog = "Non-white", var = "FearCrime") %>%
  bind_rows(results, .) -> results

corr.results[4,] <- cor.test(results$Trend[results$Demographic == "Over 50"], results$Trend[results$Demographic == "Under 50"])$estimate
corr.results[5,] <- cor.test(results$Trend[results$Demographic == "Male"], results$Trend[results$Demographic == "Female"])$estimate
corr.results[6,] <- cor.test(results$Trend[results$Demographic == "White"], results$Trend[results$Demographic == "Non-white"])$estimate


Data[[3]]  %>%
  mutate(Date = as.Date(Date)) %>%
  filter(Date > as.Date("1980-01-01")) %>%
  select(Varname, Date, Index, Poll, Demographic, n) %>%
  mutate(Demographic = ifelse(Demographic %in% c("16-24","25-49"),"Under 50",
                              ifelse(Demographic %in% c("50-64","65+"),"Over 50",
                                     Demographic))) %>%
  group_by(Date, Demographic, Varname) %>%
  summarise(Index = weighted.mean(Index, n)) -> temp

results <- analyse_data(data = temp, demog = "Over 50", var = "FearCrime")
analyse_data(data = temp, demog = "Under 50", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp, demog = "Male", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp, demog = "Female", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp[temp$Date >= as.Date('1991-01-01'),], demog = "White", var = "FearCrime") %>%
  bind_rows(results, .) -> results
analyse_data(data = temp[temp$Date >= as.Date('1991-01-01'),], demog = "Non-white", var = "FearCrime") %>%
  bind_rows(results, .) -> results

corr.results[7,] <- cor.test(results$Trend[results$Demographic == "Over 50"], results$Trend[results$Demographic == "Under 50"])$estimate
corr.results[8,] <- cor.test(results$Trend[results$Demographic == "Male"], results$Trend[results$Demographic == "Female"])$estimate
corr.results[9,] <- cor.test(results$Trend[results$Demographic == "White"], results$Trend[results$Demographic == "Non-white"])$estimate

as_tibble(corr.results) %>%
  mutate(Demog = rep(c("16-49 / 50+","male / female", "white / non-white"),3),
         var = c(rep("Fear of Crime",3),rep("Punitiveness",3),rep("Death Penalty", 3))) %>%
  pivot_wider(names_from = Demog, values_from = V1) %>%
  select(var, `male / female`, `white / non-white`, `16-49 / 50+`) %>%
  gt(rowname_col = "var") %>%
  tab_header(title = md("**Table 3: Bivariate correlations across different demographic groups**")) %>%
  cols_align(align = "center", columns = c(2:4)) %>%
  fmt_number(decimals = 2, pattern = "{x} ***") %>%
  tab_source_note(source_note = "Notes: Data from British Election Study, British Social Attitudes Survey, the Crime Survey for England and Wales, YouGov, Gallup and Ipsos. Series by age and by gender are estimated from 1981 to 2023 (n = 43). Due to small sample sizes for the non-white category in earlier surveys, those series are calculated from 1991 to 2023 (n = 33).") %>%
  tab_style(style = cell_text(style = "italic"),locations = cells_column_labels())
#gtsave(filename = "table_3.html", path = "/Users/matteo/Downloads")
rm(results, temp, corr.results, analyse_data)



#### Punitiveness by age

Data[[2]] %>%
  mutate(Date = as.Date(Date)) %>%
  filter(Date > as.Date("1980-01-01")) %>%
  select(Varname, Date, Index, Poll, Demographic, n) %>%
  mutate(Demographic = ifelse(Demographic %in% c("16-24","25-49"),"Under 50",
                              ifelse(Demographic %in% c("50-64","65+"),"Over 50",
                                     Demographic))) %>%
  group_by(Date, Demographic, Varname) %>%
  summarise(Index = weighted.mean(Index, n)) -> PunitivenessByAge

output <- DyadRatios::extract(PunitivenessByAge$Varname[PunitivenessByAge$Demographic == "Over 50"], 
                              PunitivenessByAge$Date[PunitivenessByAge$Demographic == "Over 50"], 
                              PunitivenessByAge$Index[PunitivenessByAge$Demographic == "Over 50"], 
                              unit="A", begindt=NA, enddt=NA,
                              smoothing=TRUE)
output1 <- DyadRatios::extract(PunitivenessByAge$Varname[PunitivenessByAge$Demographic == "Under 50"], 
                               PunitivenessByAge$Date[PunitivenessByAge$Demographic == "Under 50"], 
                               PunitivenessByAge$Index[PunitivenessByAge$Demographic == "Under 50"], 
                               unit="A", begindt=NA, enddt=NA,
                               smoothing=TRUE)

# Figure 5: Trends in punitiveness by age group, 1981 - 2023
plot(x = output$period, y = output$latent1, type = "l", 
     xlab = NA, ylab = "Latent trend", ylim = c(0.46,0.63),
     xlim = c(1978,2023)) +
  lines(x = output1$period, y = output1$latent1, type = "l", lty = 2)+
  text(x = c(1981, 1981), y = c(output$latent1[[1]], output1$latent1[[1]]),
       labels = c('Over 50s','Under 50s'), pos = 1)
rm(output, output1, PunitivenessByAge, csign)
